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Abstract 

We compute the improvement coefficient c sw that multiplies the Sheikholeslami-Wohlert 
term as a function of the bare gauge coupling for two flavour QCD. We discuss several 
aspects concerning simulations with improved dynamical Wilson fermions. 

1 Introduction 

The standard formulation of lattice QCD by Wilson has been used since the early days of lattice 
gauge theory. However, it is known that in this formulation the leading discretization errors 

*Talk given by K.J. at the International Symposium on Lattice Field Theory 21—27 July 1997, Edinburgh, 
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are linear in the lattice spacing a. Moreover, by testing the PCAC relation on the lattice, it 
could be demonstrated that the effects of these discretization errors are most severe and can 
influence values of physical observables strongly [[[]. 

In a series of papers [fj, |3], || it was shown that, by implementing Symanzik's improvement 
programme || for QCD on-shell and non-perturbatively, one can reach a complete cancellation 
of the 0(a) effects. The advantages of this procedure are obvious, and this conference has seen 
the improvement programme successfully at work ||. The complete improvement programme 
demands as a first step a computation of the parameter c sw that multiplies the Sheikholeslami- 
Wohlert term |7j] in the improved action. In addition, also the coefficients that enter the 
improved operators have to be determined. By now, in the quenched approximation, a number 
of these parameters are known as a function of the bare gauge coupling g j|, |8], Kl- 
in this contribution we want to initiate the computation of the improvement coefficients for 
two flavours of dynamical fermions. As a first step we will compute the coefficient c sw . As is well 
known, dynamical fermion simulations are very demanding even with today's computers and 
algorithms []10| . On the other hand, knowing the improvement coefficients will substantially 
reduce the computational cost, since one is allowed to choose larger lattice spacings. 

2 The improvement condition 

The idea of testing the lattice artefacts is to probe the PCAC relation, which should hold, up 
to 0(a) corrections: 

d^x) = 2mP a (x) + 0(a) , (1) 

where A a (x) denotes the isovector axial current and P a (x) the corresponding density. The 
quark mass that appears in eq. (|XJ) is a bare current quark mass at scale 1/a. The important 
point to notice here is that the PCAC relation is an operator identity that can be inserted into 
arbitrary correlation functions. 

One can make use of this fact to improve the theory: one tests the PCAC relation in 
different correlation functions and demands to obtain always the same value of m. Using the 
Schrodinger functional, it was demonstrated in ref. [||] how this strategy can be efficiently 
implemented to determine c sw . Here we follow ref. |Q closely and impose exactly the same 
improvement condition: 

aAM = 0.000277 (2) 

at L/a = 8, with AM as introduced in ref. [|J. 

The improvement condition eq. (§) can in principle be imposed at any (not too large) value 
of M, where M is a specific definition of the current quark mass [[| derived from eq. (1). In 
order to guarantee a smooth behaviour of c sw (g ) one should, however, make a definite choice, 
where a natural value is M = 0. Since AM turns out to be a very weak function of the quark 
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mass M 0, |TT|, one may also compute AM for \aM\ <C 1. In particular, for the values of aM 



chosen, here, the error introduced is negligible compared with the statistical one. 

3 The simulations 

The numerical simulations are performed on 16 x 8 3 lattices, with boundary conditions as 
detailed in ||. We use the Hybrid Monte Carlo (HMC) algorithm with the Sexton- Weingarten 
scheme to integrate the classical equations of motion f!2]| . Our implementation of the HMC 



algorithm is described in detail in ref. [13]. All simulations are performed on the massively 
parallel Alenia Quadrics (APE) computers. On the two versions of these machines that we have 
used with 256 and 512 nodes, we ran 32 and 64 independent simulations in parallel. Combined 
with a jack-knife method, this allows for a realistic error estimate on our observables. 

We have run simulations at eight values of (3 = 6/(?q in the range 5.2 < (3 < 12.0. Each 
simulation has at least 1280 molecular dynamics trajectories and typically 2500. Keeping the 
trajectory length fixed to 1, we reach typical acceptance rates of 95%. Despite the relatively 
large acceptance rates, we noticed that sometimes a system can get stuck and does not accept 
a number of (larger than, say, 10) trajectories. The problem is easily overcome by performing 
every n number of trajectories one with a much smaller step size. Of course, in order to 
be able to show that one generates the correct distribution, the value of n has to be chosen 
independently of the Monte Carlo history. We simply kept n fixed in each simulation. 

We applied the improvement condition eq. (0) in the small quark mass region, \aM\ < 0.01 .Q 
With Schrodinger functional boundary conditions, simulations at such small quark masses are 
unproblematic, since the massless Dirac operator of the Schrodinger functional with time extent 
T has a lowest eigenvalue with mag nitude A min = const. /T + 0(g 2 ) + O(M). 

Owing to the 0(g 2 ) terms in A m ; n , the simulations slow down, when f3 is decreased. In detail 
the reason for this is threefold. First, going to smaller values of f3 we have to decrease the step 
size dt from, as an example, dt = 0.066 at (3 — 7.4 to dt = 0.027 at (3 — 5.4. Second, the 
condition number k of the preconditioned fermion matrix Q 2 (see e.g. ref. [[HI for a definition 



of Q) increases with decreasing (3, as can be seen in fig. 1. The increasing values of k result 
in a growing number of conjugate gradient (CG) iterations when going to smaller (3. Third, 
we find an increase of the autocorrelation time r with decreasing (3 for observables such as the 
lowest eigenvalue of Q 2 or quark correlation functions at a given distance. Fortunately, it turns 
out that the autocorrelation time for AM is small, r m 2-4, and shows only a weak dependence 
on (3. 

As a rule, we find that the performance of the simulation algorithm does not significantly 
depend on the value of c sw . The only exception are the autocorrelation times r for which there 
are indications that they are particularly large when both c sw and (3 are small. 

1 There is one data point where aM = 0.023. 
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Figure 1: The condition number k of the preconditioned fermion matrix Q 2 . Different values 
of k at the same bare gauge coupling correspond to different values of c sw . 

4 Results 

We determined AM at fixed (3 for various values of c sw . From the -linear- dependence of AM 
on c sw we can then extract the slope s = dAM/dc sw . We found in practice that the slope is well 
described by a linear function of g%. In order to extract the desired improvement coefficients 
Cg™ pr (#o) at the eight values of g l ,i — 1, ...,8, where the simulations are performed, we fit all 
our data for AM to the form 



where 



a AM = s(g ) ■ (c sw - C^*)) = 0.000277, 



s(g ) = -0.015 • (I + Sl gl) 



(3) 



(4) 



and Si as well as c l ™ pT (gl) are fit parameters. The results for Cg™ pr are displayed as the full 
symbols in fig. 2. The solid line is a representation of these data, given by 



I - 0.454^ _ 0.175^ + 0.012^ + 0.045^ 



0.720^, 



(5) 



As already mentioned, for small values of (5 the simulations become very costly, and we were 
not able to perform simulations at (5 — 5.2 and small quark masses. We therefore switched to 
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Figure 2: The improvement coefficient Cg^r CIS 9b function of the bare gauge coupling g^. The 



solid line represents eq. (J5|). The dotted line is the 1-loop result || [14| and the dashed line is 
the result in the quenched approximation H. 



the following strategy: we take the parameterization eq. @ to extrapolate a little bit further 
in p,to/3 = 5.2. At the value of c sw determined in this way, we then select a large quark mass, 
aM = 0.1 and try to verify that improvement is at work. Indeed, we find for f3 = 5.2 and 
c sw = 2.02 that aAM = —0.0006(9). This indicates that our final result eq. (|5|) can safely be 
used for (3 > 5.2. Preliminary studies of the hadron spectrum in the improved theory suggest 
that (3 > 5.2 yields the range of lattice spacings that is of interest to computations of hadronic 
properties [P"5| . 

We want to emphasize that although, with our values of c sw , the O(a) terms are cancelled, 
0(a 2 ) effects remain and are not negligible for (3 ~ 5.2, as will be discussed elsewhere [11 . 



This work is part of the ALPHA collaboration research programme. We thank DESY for 
allocating computer time to this project. 
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